Dynamical instability in a trimeric chain of interacting Bose-Einstein condensates 
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We analyze thoroughly the mean-field dynamics of a linear chain of three coupled Bose-Einstein con- 
densates, where both the tunneling and the central-well relative depth are adjustable parameters. Owing 
to its nonintegrability, entailing a complex dynamics with chaos occurrence, this system is a paradigm for 
longer arrays whose simplicity still allows a thorough analytical study. We identify the set of dynamics 
fixed points, along with the associated proper modes, and establish their stability character depending on 
the significant parameters. As an example of the remarkable operational value of our analysis, we point 
out some macroscopic effects that seem viable to experiments. 



5h 



c 

o 
o 



> 

m 
in 
o 

O 

-i— > 



T3 

c 

o 
o 



X 



PACS numbers: 03.75.Fi, 05.45.-a, 03.65.Sq 

Since the first realizations of Bose-Einstein conden- 
sation in atomic gases, great efforts have been aimed at 
improving the control of geometrical arrangements of the 
condensate and, in particular, at realizing its fragmenta- 
tion in many interacting components. The design of in- 
creasingly efficient trapping schemes [j, |4J has shown 
this program to be a fairly realistic perspective. Indeed ex- 
periments involving Bose-Einstein condensates (BEC) dis- 
tributed within optical traps were successfully conducted 
that provided quite large ID |5] and 2D [6] periodic arrays 
interacting via tunneling. Based on magnetic trapping I7J] 
the simplest case, consisting of two coupled condensates 
(dimer), was realized as well. In parallel, the rich scenario 
of phenomena observed in BECs arrays (nonlinear oscilla- 
tions Opl, self-trapping |0], supercurrents flfjll . and solitons 
ITTII12II ') raised a number of questions on their time evolu- 
tion. The possibility to detect and study, both experimen- 
tally and theoretically, new macroscopic dynamical phe- 
nomena/effects in BEC arrays (thus getting a deeper insight 
of stability properties and, operatively, an increased control 
of systems) has prompted an intense ongoing work. 

In this perspective, the three-site array (trimer) deserves 
a special attention. According to the standard mean-field 
treatment, which is fairly satisfactory when the average 
well populations are large fl3ll . the dimer dynamics is in- 
tegrable 1 14]. The latter is described by two macroscopic 
complex variables Zi = \zi \ exp(i#j), accounting for the 
condensates' state (phase and population |^| 2 ), and ex- 
hibits two constants of motion, namely the total boson 
number and the energy. The apparently harmless addi- 
tion of a further coupled condensate is sufficient to make 
the system nonintegrable thus causing, in the presence of 
nonlinear BEC self-interactions, strong instabilities in ex- 
tended regions of the phase space 1 15]. That is, while keep- 
ing simple enough to be viable to a thorough analytical 
study, the system displays a whole new class of behav- 
iors which, though often overlooked, are typical of longer 
arrays. Moreover, the recent achievements in the experi- 



mental field - and especially the control promised by mi- 
cro traps 1 4] - suggest the realization of the trimer to be at 
hand. For all these reasons we feel that the trimer deserves 
a systematic analysis, whose key results we discuss below. 

We consider an asymmetric open trimer (AOT) made of 
three coupled BECs arranged into a row, where both the 
interwell tunneling T and the central-well relative depth w 
are adjustable parameters. Compared with the symmetric 
case iflal (trimer on a closed chain of equal-depth wells), 
the interplay of such parameters both entails a deeper con- 
trol on the dynamics and favours the approach to experi- 
mental situations. After recognizing the location of fixed 
points in the phase space and the associated proper modes, 
we focus mostly on establishing via standard methods their 
stability character on varying r = T/UN and v = w/UN, 
where U embodies the interatomic scattering and N is 
the total boson number. Also, we show how macroscopic 
(i.e. interesting experimentally) dynamical effects can be 
primed by selecting suitable critical values of r and v. 

The essential physics of the AOT is aptly described fl3l 
Il7ll by the Bose-Hubbard Model, which represents a gas 
of identical bosons hopping across an ambient lattice. The 
relevant Hamiltonian, ensuing from the boson field theory 
through the space-mode approximation of field operators 
fill, reads H = E| =1 {Un\ - vn k ) - wn 2 - T[a^(a 1 + 
a 3 )+h.c.]/2, where v is the depth of wells j = 1, 3 (w, T, 
U are defined above), = a\ai counts the bosons at site i, 
and the destruction (creation) operators eij (af) obey com- 
mutators [a,i, a~l] = 5ih- If the well populations are not so 
small that a purely quantum treatment is in order, the sys- 
tem dynamics can be described by three variables Zi within 
the mean-field picture. This is recovered via the coherent- 
state variational procedure of Ref. 17, where the semiclas- 
sical Hamiltonian Ti = T 1 \ =l {U\zk\' i: — v\z k \ 2 )—w\z2\ 2 — 
T/2[z2 (zi + z 3 ) + c.c.] is derived from the system effec- 
tive action S(t) ensuing in turn from the macroscopic trial 
state = exp (iS/H) Hi \zi), written in terms of the 



Glauber's coherent states, a, 



In this frame- 
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work the symbol n$ is used to denote the expectation value 
(\P|ni|\l/) = \zi\ 2 . Equipped with the Poisson brackets 
{z* ,z h } = i/h5j h ,TC yields the dynamical equations, 

( ihzj = (2U\z j \ 2 - v) Zj - § z 2 m 
\ ihz 2 = (2U\z 2 \ 2 - v - w) z 2 - ^ ( Zl + z 3 ) ' W 

with j = 1, 3 in the first equation. Eqs. Q incorporate the 
conservation of total boson number N = T,j =1 \zj\ 2 , since 
{N, TC} = 0, and provide, by conjugation, the equations 
for {zl}. Further, they propagate in time the initial condi- 
tion Zi = z 3 , thus revealing an integrable subregime char- 
acterized by the first integrals H and N = 2\zi\ 2 + \z 2 \ 2 : 
We refer to such a regime as the dinieric regime (DR), as 
opposed to the nondimeric regime (NR) where Z\ ^ z 3 . 

The distinctive features characterizing a given Hamilto- 
nian dynamics are deduced by exploring its phase-space 
structure. This is attained in the first place by working out 
the location and the local character of its fixed points. For 
Eqs. 0, the latter are issued from 

f = (2U\z ] \ 2 - f i)z J -^z 2 j = l,3 

\ 0= (2U\z 2 \ 2 — fi — w) z 2 -^( Zl + z 3 ) (1) 

where /i = v + x and X is a Lagrange multiplier selecting 
the conserved value of N. Every solution (r]i,7] 2 ,rj 3 ) of 
Eqs. (0 naturally provides a periodic solution of Eqs. (0 
of the form Zj(t) = Tjj exp[ixt/H]. The global phase 
symmetry Zj i— > Zj exp[i<3?] ($ G [0, 2ir]) allows one to 
replace z t with Xi £ R thus reducing system to 

f 0={2x 2 -m)x J -\x 2 J = 1,3 

X § = {2x\-m-v)x 2 -\ (xi + x 3 ) 

where m = fx/UN, and x* = (xi,x 2 ,x 3 ) is such that 
Zj = \/Nxjexp[i§>] and S i a:? = Sjni/iV = l. Notice that 
the structure of system © entails that the solution where 
either one or both of the peripheral sites are depleted is 
trivial: x X 3 = =>■ x = 0. The situation where the central 
well is depleted, henceforth referred to as central-depleted 
well (CDW) yields x x = -x 3 = ±y/m/2 = ±l/y/2. 
Note that this represents the trimeric counterpart of dimer 
7r-states 1 8]. In general, except for the simple CDW case, 
solving system (|3j for r ^ is found to be equivalent to 
finding the real roots of the quartic polynomial in a [20] 

V TV {a) = a 4 + b T „a 3 + c T „a + d = 0, (4) 

Once a solution is found in terms of a, the relevant config- 
uration x is recovered as x 2 = aR(a) /^l + a 2 , 

x lt3 = (X x ±X 3 )/>/2 (NR), x h3 = X h3 (DR) (5) 

where X 3 = (1 - X\ - xf)i, X 1 = R(a)/^l + a 2 
and R(a) = (1 - ra/\/8)^ in the NR, while X 1 = 
R(a)/ V / 2(1 + a 2 ), R(a) = 1 in the DR. 

The study of Eq. © allows one to construct the upper- 
panel diagrams of Figs. ^ an d where the same shade 
of gray characterizes regions of the tv plane featuring 
the same number of fixed points. Also, it represents an 



essential element in building the stability diagrams [dis- 
played in the lower panel of Fig. ^ (Fig. |2jl for the DR 
(NR)] which are the central result of this letter. Opera- 
tionally, such diagrams, depicted in the 9t plane (where 
9 = arctan(a) G [— vr/2, vr/2]), allow to determine both 
the location and the stability character of the fixed points. 
For a given parameter pair (r, v), the latter are obtained 
by intersecting the straight line r = f and the curve 
t = Tp(6), where r = T u (6) is implicitly defined by 
equation V tu {ol) = 0. The fixed-point components a;j's 
are obtained from the relevant 9 via Eqs. The lower 
panel of Fig. ^l|2j displays the graph of r = t v {9) for 
five (four) choices of v, which we discuss below. We re- 
mark that such diagrams embody the whole information 
on trimer dynamics. There, the fixed points are identified 
as maxima, minima and (stable/unstable) saddles based on 
the shade of gray of the region they lie in. 

Fixed-points local character. This is identified ana- 
lytically through a standard study of the signature of the 
quadratic form associated to a given Hamiltonian, along 
the lines of Ref. Il6[ Locally, to the second order in the 
displacements £ 3 - = \[N (qj + ipj) = z 3 — C, 3 around the 
fixed point Q = yNxj, one gets 7i(Zj) — xN =TC((j) + 
U N 2 t h(q,p)/2, with(i t = (q 1 ,q 2 ,q 3 ),p t = (p 1 ,p 2 ,p 3 ) 
and /i(q, p) = q t Qq + p*Pp. Remarkably, the six-by-six 
matrix relevant to the quadratic form /i(q, p) in the vicin- 
ity of a fixed point, separates into two three-by-three sym- 
metric submatrices P, Q depending on the Xj's, t and v. 
Also, when the displacements' constraint x*q = (issued 
from Ej | z 3 \ 2 = N) is accounted for, the rank of Q fur- 
ther reduces to two. Likewise only two eigenvalues of P 
affect the fixed-point character, the third one being iden- 
tically zero. Indeed it is easy to check that the relevant 
eigenvector is p oc x. The non-zero eigenvalues of P have 
opposite signs in the CDW configuration (which is there- 
fore a saddle), and the signs oix\/x 2 and x 3 / x 2 otherwise. 
Hence the character of a fixed point depends on the signs 
of the two pairs of significant eigenvalues of P and Q. 

Linear stability analysis. The linear stability 

HI of 

the fixed point can be determined by studying the evo- 
lution of the displacements q and p. This is governed 
by the differential equations q = {q,7i} = aPp and 
p = { P) H} = -o-Qq (with a = UNt/2K), through 
the Poisson brackets {q 3 ,Pk} = ^jfe/(2 NK) stemming 
from the original ones {z* , = iSj^/h. Such qp equa- 
tions provide the dynamical six-by-six matrix M 6 with di- 
agonal (three -by-three) blocks vanishing, and off-diagonal 
blocks that coincide with matrices P and — Q 12211 . 

Based on the standard criteria, a fixed point is (linearly) 
stable when none of the six complex eigenvalues {\j}® =1 
inherent in M 6 features a positive real part. Otherwise the 
fixed point is unstable and a chaotic behavior may arise. 
The Ai's are conveniently obtained as the square roots of 
the eigenvalues {A^}| =1 of —PQ. Indeed, due to the block 
off-diagonal form of M 6 and to the symmetry of P and 
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FIG. 1: Stability diagram for the DR. Upper panel: number of 
fixed points. Lower panel: character of the fixed points. The 
color keys of both panels, sharing the same r axis, are the same 
as in the GR case (see caption of Fig. [2}. 



Q, det(A-M 6 ) = det(A 2 + PQ). As expected, this im- 
plies that maxima and minima are stable points. Indeed for 
such configurations the eigenvalues of P and Q have the 
same signs (non-positive for the maxima and non-negative 
for the minima), yielding therefore non-positive Aj's and 
purely imaginary Aj's. Further, the fact that one of the 
eigenvalues of P is zero implies that one of the A^ 's, and 
hence two of the Xj 's, are zero as well. The four remaining 
A/s are the roots of a biquadratic polynomial of the form 
A 4 — s A 2 +p. A fixed point is therefore stable if the condi- 
tions s 2 — 4p < 0, p > and s < are simultaneously met. 

Discussion. The character/stability of fixed points are 
fully described by the lower panel of Fig. \l\ (Fig. - 
referred to the DR (NR)- through the graphs of function 
T u (9) for five (four) significant choices of v, distinguished 
by different dashing patterns. In the figure upper panels a 
straight horizontal line featuring the same dashing style as 
the corresponding r v {9) curve allows one to read the rel- 
evant value of v. As noted above, each pair (v, f ) selects 
the set {9 r : f = T p (9 r ), r < 4} giving the fixed-point 
components. Depending on the shade of gray filling the re- 
gion it lies in, a fixed point is either a minimum (white), a 
stable/unstable saddle point (light/medium gray) or a max- 
imum (dark gray). In the DR (Fig. each curve r„(0) 
has two branches featuring asymptotes at 9 = ±7r / 4. For 
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FIG. 2: Stability diagram for the NR. Upper panel: number of 
fixed points. The darker is a region, the larger is the relevant 
number of solutions, also displayed. The lines at a fixed value 
of v correspond to the curves r v iff) appearing in the lower panel 
with the same dashing style. Lower panel: character of the fixed 
points. Dark gray: maxima; Medium gray: unstable saddles; 
Light Gray: stable saddles (only in the NR); White: minima (only 
in the DR); Diagonal lines: forbidden regions (only in the NR). 
Both panels share the same r axis. 



— 1 < v < 2, a third bell-shaped branch [see T u (9) for 
v = 0.2, 0.4] crops up. Hence the DR always has two 
fixed points, of which one is always a minimum, whereas 
the other is a maximum (unstable saddle) for large (small) 
r's. For v = 1/2 the unbounded branches actually col- 
lapse on their asymptotes, thus providing r-independent 
solutions |23J|. In the NR, Fig. |2] t„{9) features one or 
two bell-shaped branches, depending on whether \u\ > 2 
or \u\ < 2, respectively. Notice that in some cases entire 
portions of such branches, namely the ones lying within 
the patterned regions, must be discarded since, despite they 
are real, the relevant roots of polynomial © yield complex 
solutions of system 0. The dotted curves appearing in 
the upper panel of Fig. are actually the counterparts of 
the dotted curves delimiting such forbidden regions, and 
significantly modify the simpler picture ensuing from the 
mere study of the real roots of polynomial @ (solid black 
lines). For r < 2[2/3 1 ^ 2 — l] 1 / 2 and v ranging around 
v = 2(3 1 ' 2 — 2) a four-solution lobe is found. As observed 
above, maxima and minima are predicted to be stable fixed 
points. Noticeably, the stability analysis evidences that in 
the NR a significant fraction of the saddle points is dynam- 
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FIG. 3: Sets S and M of Poincare sections of trajectories based 
close to Co for two choices of significant parameters r, v (see the 
text for further detail). Sections are cut at £ 2 = —3/4 within the 
reduced phase space (£i , £2 , <t>x , ^2). [with = ( n 3 — n i +^2 )/iV , 
£2 = (713 — rci — n2)/An ensuing from AT = (see also Ref. Uffl . 



ically stable (light gray), while in the DR the saddle points 
are always unstable (medium gray). Also, while Fig. ^ 
shows how the minimum always belongs to the DR, the 
absolute maximum is found either in the DR or in the NR 
depending on the values of v, r. As to the CDW saddle, it 
is unstable in the region {v< — 1,t< — z^}U{|^| < 
1, T > [(1 + v) 3 /(2 - 2 u)]i} and stable elsewhere. 

As already mentioned, the versatility of micro traps |@1 
promises a full control on the features of the trapping po- 
tentials, and hence on the parameters r and v. It is worth 
noting that present experimental setups o,|6J] already allow 
to tune the ratio T/U so that, according to the estimates in 
Ref.E r spans the range displayed in Figs. U and |2] if a 
macroscopic N = 10 3 -f- 10 5 is considered 12411 - 

We conclude by emphasizing two interesting features of 
AOT, among the many suggested by our analysis, which 
seem to be more readily apt to experimental tests. As to 
the just mentioned CDW saddle we recall that its form 
(n 2 = 0, counter-phased side wells with n\ = n 3 ) does 
not depend on parameters r and v. Hence a system pre- 
pared in such a configuration with a weak (r <C 1) inter- 
well tunneling, would remain in such apparently unnatural 
state also when tunneling is enhanced (r 3> 0), provided 
the central well is sufficiently deep (y > 1) to ensure the 
dynamical stability of the fixed point. Notice that preparing 
two separate counterphased condensates is by no means an 
easy task. In this respect the above described phenomenol- 
ogy suggests an operational method to check if the system 
is satisfactorily close to the desired configuration. Further, 
we remark that, owing to the relatively simple algebraic 
structure of the dimeric equations, there exists a simple, 
actually linear, functional relation among r and v ensur- 



ing that a fixed point is characterized by the same coordi- 
nates x.i's independent of the parameters [25]. This allows 
one to evaluate with little effort different sets of parameters 
(r, u) for which the same fixed point belongs to regions 
having a different stability character. Then, a relatively 
simple tuning of r and v leads to change the fixed-point 
character from stable to unstable, thus inducing chaos on- 
set. Such an effect is manifest in Fig. where two sets 
(labeled by M and S) of Poincare sections, both issued 
from initial conditions quite close to the configuration C 
(m =n 3 = N/8, n 2 = 3/ 'W, fa = # 2 - tf x = vr, fa = 
"9z-$2 = -vr) areplottedfor (r M , z^m) = (7/[4^6], 2/3) 
and (t s ,u s ) = (13/[8\/6], 17/24), so that the corre- 
sponding points (t m s ,# = — vr/3) lie just in the max- 
ima and the saddle regions, respectively. Notice indeed that 
the former parameter choice yields the regular section ex- 
pected for a stable point. Conversely, when the fixed point 
is a saddle, trajectories based very close to it invade a large 
phase-space region, as testified by Fig. where the rele- 
vant Poincare sections are undistinguishably merged into a 
fuzzy cloud of points. The analytical study exposed in the 
present letter is widely confirmed by numerical simulations 
where the chaoticity of the unstable fixed points is mani- 
fest. A thorough discussion of the AOT phenomenology 
and of further macroscopic effects disclosed by the stabil- 
ity diagrams will be presented in a later paper [24]. 
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